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ABSTRACT 



Aims. The methodology for modeling grain-surface chemistry has been greatly improved by taking into account the grain size and 
fluctuation effects. However, the reaction rate coefficients currently used in all practical models of gas-grain chemistry are inaccurate 
by a significant amount. We provide expressions for these crucial rate coefficients that are both accurate and easy to incorporate into 
gas-grain models. 

Methods. We use exact results obtained in earlier work, where the reaction rate coefficient was defined by a first-passage problem, 
which was solved using random walk theory. 

Results. The approximate reaction rate coefficient presented here is easy to include in all models of interstellar gas-grain chemistry. 
In contrast to the commonly used expression, the results that it provides are in perfect agreement with detailed kinetic Monte Carlo 
simulations. We also show the rate coefficient for reactions involving multiple species. 
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1. Introduction 

The chemistry in interstellar clouds of gas and dust is very com- 
plex, consisting of reactions takin g place in the gas phase as well 
as on the surfaces of dust grains (Herbsn ll995h. Of the latter re- 
action s, the formation of molecular hydrogen dHollen bach et al. 
1 197 lh is of parti cular importance, but m any more processes have 
been identified (Hasegawa et al. 1992). Modeling of these com- 
plex gas-grain chemical networks is important to the understand- 
ing of the observed chemical complexity in interstellar clouds 
and the effects of those molecules on gravitational collapse and 
star formation dHerbstll995t lHartguist & Williamsll995t(Tielensl 
I2005h . 

While for gas-phase reactions, rate equations are widely 
used and fully appropriate, for grain surface chemistry they 
were found to be ina ccurate for a wide range of astroph ysically 
relevant condition s dTielensI 1 1 9951: ICaselli et al.lll998k Biham 
& Lipshtat 120021: iBiham et al] 120051: lLohmar & Krugl 120061: 
, iLederhend ler & Biham 2008). To address this proble m, modified 
rate equations were proposed ICaselli et alj 19981: Stantcheva 
' et al. 120011 : iGarrodl l2008h . Thev were found to provide im- 
proved results, but involve ad-hoc and uncontrolled approxima- 
tions. The master equation, describing the time evolution in the 
probability distribution of the number of reactants on the grain, 
provides a com plete and accurate description of the re actions on 
grain surfaces (Biham et al . 2001 1 iGreen et al.lb OQl: Biham & 
Lipshtat l2"002h . 

The surface chemistry was also studied using stochastic sim- 
ulations dTielens & HageiJl982tlCharnlevl2 001). Kinetic Monte 
Carlo simulations (KMC ) have been developed for different su r- 
face morphologies (IChang et alJl2005l : ICuppen & Herbstll2005h . 

An important advantage of approaches based on differential 
equations is that they can easily be coupled to the rate equations 
of gas-phase chemistry. The difficulty with the master equation 



is that the number of equations quickly proliferates as the num- 
ber of reactive sp ecies on the grains increases, thus approxima- 
tions are needed dStantcheva et al.ll2002t[Stantcheva & Herbst 
2003). From the master equation, it is possible to construct a 
set of moment equatio ns that account correctly for the reaction 
rates o n grain surfaces (Lipshtat & Biham 2003; Barzel & Biham 
l2007allbh . 

Here we consider the form of the rate coefficients that en- 
ter the models. Apart from KMC simulations, all approaches to 
grain-surface chemistry feature the sweeping rate A as a crucial 
parameter. This rate coefficient governs the reaction rate on the 
surface, where there is competition between diffusion-mediated 
reaction and desorption. Without a precise definition of A, the 
expression A = a/S is widely used, where a is the hopping 
rate between adjacent adsorption sites and S is the number of 
these sites on the grain. This expression does not account for the 
competition between reaction and desorption. It also does not 
account for the peculiarities of two-dimensional diffusion. 

Lohmar & Krug (2006) introduced a proper definition of the 
rate coefficient A in terms of the two-particle encounter prob- 
ability. These authors evaluated it exactly for a single-species 
reaction on a closed two-dimensi onal surface where the adsorp - 
tion sites form a regular lattice (Loh mar & Krugll2006l 12009). 
The result was compared to the conventional approximation for 
A, which was found to deviate significantly. However, the exact 
expression contains a large sum of terms, which is costly to eval- 
uate in practice. 

In this Letter, we present an expression for the reaction rate 
coefficient that is both accurate and easy to evaluate. We com- 
pare the reaction rates obtained using this expression to those 
obtained from KMC simulations and find perfect agreement. We 
also generalize the expression to the case of reactions involving 
multiple species. 
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2. Review of earlier results 



3. Practical approximation 



There is no need here to elaborate on the various analytical ap- 
proaches to grain-surface reaction kinetics, which we cited be- 
fore. Throughout, we consider a system of S adsorption sites 
arranged in a two-dimensional quadratic square lattice with pe- 
riodic boundaries, onto which atoms impinge at a certain homo- 
geneous flux F — fS . They perform nearest-neighbor hops with 
(undirected) rate a and desorb at a rate W <s a. The sweeping 
rate A then appears in all approaches as the pre-factor in the re- 
combination term and the corresponding production rate. For ex- 
ample, the rate equation for the mean number of particles on the 
grain, (N), reads d(N)/dt = F - W(N) - 2A(N) 2 , with a produc- 
tion rate A(N) 2 of molecules. In the master equation and the mo- 
ment equations, the latter rate becomes instead R = A((N 2 )-(N)) 
dBiham & Lipshtatll2002l; iLipshtat & Bihamll2003l e.g.), where 
the moments (N k ) = 'E l ^ =0 N k P(N), and P(N) is the probability 
of having N particles on the grain at a given time. Any adsorbed 
atom may end up either reacting (with rate A) or desorbing (with 
rate W). Therefore, the probability for an atom to end up react- 
ing is giv en by p — A I (A + W). With this, we can define A by 
means of dLohmar & K rug 20061) 



A 
W 
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(1) 



where p is the encounter probability that two specific particles 
present on the grain meet before one of them desorbs. If atoms 
only react with a certain probability upon encounter, p is simply 
multiplied by this probability. In Sect. [5] below, we derive the 
corresponding relation in the context of multi-species reactions, 
and also show how the encounter probability generalizes to that 
case. 

The encounter probabilit y has been studied in great detail 
dLohmar & Kru g 2006, 2009). Generally, it features two regimes: 
the 'small-grain' regime, where SW/a <K 1 and encounter is 
almost sure, 1 - p <k 1, and the 'large-grain' regime, where 
SW/a s> 1 and encounters are rare, p <K 1. The best model 
at hand is given by a discrete-time random walk analysis. Well- 
known results then yield p in t erms of a single sum, and the 
asymptotic behavior is given bv dLohmar & Krug||2009h 
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(with a numerical constant c = 1.8456047840 . 
grain regime, while in the large-grain regime, 
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.) in the small- 
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Using Eq. ([]]) thus yields the small-grain asymptotic 
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and for large grains, 
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S ln(8a/W) 



(5) 



This behavior is largely independent of the precise model used 
to obt ain p, the latter only affe cting the pre-factor inside the loga- 
rithm dLohmar & Krugll2009h . The logarithmic factor in Eqs. (|4]i 
and (O reflects the return probability of a two-dimensional ran- 
dom walker ("back diffusion") and as suc h, is also t ypical of 
other two-dimensional diffusion problems dKrudl2003h . 



The exact expressions for p are typically complicated, but the 
asymptotic results for the sweeping rate of interest are fairly sim- 
ple. Moreover, its functional form is the same in both regimes, 
only differing in the argument to the logarithm, which naturally 
suggests an interpolation formula. 

We first considered to simply use the minimal argument to 
the logarithm, making the dependence continuous and produc- 
ing the correct limiting behavior in both regimes. While this 
appears to be natural and is very easy to implement, it also in- 
troduces a cusp to the dependence, and it is an arbitrary choice: 
Any function of cS and 8a /W that approaches the smaller argu- 
ment in the corresponding limit could be used as well. We thus 
tested [(cS)~ n + ($a/W)~"]~ 1/n instead. Evidently, this has the 
proper limits for n > 0, and the exponent can be used to tune 
the behavior in-between; the minimum function corresponds to 
n — > oo. For several values of a/W » 1, we "visually" tuned 
the exponent to yield the smallest relative deviation from the 
exact values, and ended up at values n = 1.07 .. . 1.08, which 
decrease the maximum deviation with respect to the exact result 
from between 3-7% (for the min function) to ~ l/30th of this 
value. Numerically, this however is an extraordinary complica- 
tion. Being so close to unity, n = 1 is a feasible choice instead, 
and this still decreases the maximum deviation by a factor of 
about 1/6 in the examples shown below. Consequently, we pro- 
pose the approximation 



A ~ — • 7T 

S 
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1 



l/(cS) + W/(Sa) 



(6) 



where c = 1.8456047840 This is the main result of this 

Letter, and we advocate its use in all models of grain-surface 
chemistry. 

3.1. Comparison 

In Fig.Q~l we present the exact sweeping rate together with the ap- 
proximate expression in Eq. (0, both normalized by the conven- 
tional approximation a/S . The approximation proposed above is 
in perfect agreement with the exact result. Using thermal activa- 
tion of rat es with activat ion energies for amorphous carbon as 
given bv lKatz et all d 19991) . W/a = 10~ 3 occurs at T ~ 21.34K, 
and the other values at 16.00 K and 10.67 K, respectively. The 
discrepancy with respect to the conventional approximation is 
found to be most pronounced for large grains. We note that even 
for the smallest grains, where the correction is least severe, the 
factor is still well below unity. 



(4) 3.2. Lattice-type effects 



In the above, we considered the case of a quadratic square lat- 
tice. The restriction to the quadratic case, in which the lattice 
extends to the same lengths in both dimensions, is not easily 
removed; indeed, rectangular lattices ex hibit a more complex be - 
havior than desired for a simple model dLohmar & Krugl l2009). 
Different lattice types can be examined rather easily however, 
and we also considered the triangular lattice. For the exact ex- 
pressions, the effect on the encounter probability p as well as 
the sweeping rate A has been shown to be minute throughout 
dLohmar & Rrugj|2006l l2009h . In the approximation presented 
herein, the triangular lattice simply amounts to a change in the 
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Fig. 1. Ratio between the correct sweeping rate A and the 
tional approximation a/S , as obtained from the exact expression 
for A (solid lines) and from the approximate form (circles) given 
by Eq. ©, for three different values of W/a «c 1. 



pre-factors in Eq. (O, namely 
2n 



a 



In 
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l/(c5) + W/(12a) 
where now c = 2.3472914383 .... 

3.3. Langmuir-Hinshelwood rejection 

Expressions for the encounter probability given so far include 
coinciding initial positions of both particles. The question of 
whether Langmuir-Hinshelwood rejection should be accounted 
for on the 'global' level of the surface-chemistry model is be- 
yond the scope of this Letter; most likely, it is relevant only for 
H and D atoms. Irrespective of this, we are interested here in the 
most accurate numerical value for the reaction rate coefficient 
as opposed to a consistent theoretical framework (e.g., one may 
even use rate equations in the end). Consequently, one might 
want to use 

Sp-l 



1 



<P 



(8) 



for the encounte r probability, which do es not permit coinciding 
initial positions (Lohma r & Krudl2009h . The implied change in 
the sweeping rate is given by 



A = 



Wp 



= A - 



W 



s(i- P ) 



(9) 



We evaluate this for both regimes separately. For small grains, 
W/(l - p) ~ A, and thence, A m A(l - 1/5). For large grains, 
W ~ A/p, such that 



A « A 1 



1 



Sp(i-P)) 

Here, p <k 1, so p(l - p) m p, and from Eq. (f3j), 

W\n&a/W\ 
A as A 1 ' 



(10) 



(11) 



Since both S as well as a/W can safely be assumed to be larger 
than a few hundred, the relative corrections to A are below 1%, 
and this is as good as can be achieved with the interpolation in A. 
Consequently, accounting for Langmuir-Hinshelwood rejection 
in the sweeping rate in practice is unnecessary. 



4. Kinetic Monte Carlo simulations 

We now show that using the approximation in Eq. ((6} in analyti- 
cal approaches is in excellent agreement with KMC simulations. 
In our context, we restrict ourselves to a proof-of-principle com- 
parison, hence we treat only the case of single-species recombi- 
nation. We also use the master equation approach as the most pre- 
cise model, to avoid discrepancies that might arise from factors 
unrelated to our approximations. In practice, one would use an 
appro ximation such as the moment equations (Lipshtat & Biham 
2003), which have been shown to produce accurate results for a 
broad range of conditions, and to be applicable to large-scale 
astrochemistry simulations dBarzel & Biham 200720). 

The single-species master equation admits a stationary solu- 
tion dGreen et al.l200ll : lBiham & Lipshtatl2002h . The efficiency tj 
of the process is then given by the production rate R (molecules 
produced on the grain per unit time) normalized by the influx, 
and one obtains 



R 



i] = 



'W/A+l 



(2V2T/A) 



F/2 Iw^il^lWjA)' 



(12) 



where I^iz) is the modified Bessel function of the first kind. 

Ou r KMC simulations use a standard algorithm (described, 
e.g., in lChang et al.ll2005h . which has been suitably optimized 
for the peculiarities of our system, and the code has been tested 
extensively. Atoms impinging onto occupied sites are rejected 
by the Langmuir-Hinshelwood mechanism, while we do not ac- 
count for this process in the analytical treatment — this illus- 
trates that the framework presented here is fully appropriate for 
moderate coverages (i.e., for the high-temperature drop in effi- 
ciency). We also wait for a long enough time to establish steady- 
state conditions, and follow the system for 10 8 impingements. 

Surface parameters are c hosen as found e xperimentally for 
amorphous carbon surfaces dKatz et al.lll999h . with thermal ac- 
tivation energies £V/^b ~ 698 K, E a /k^ m 511 K, and an at- 
tempt frequency v = 10 12 s _1 . With experimentally determined 
site densities and assuming standard gas phase conditions (hy- 
drogen atom density of «h = 10 cirr 3 at a temperature of 100 K), 
one o btains an H influx / = 7 . 3 x 1 0~ 9 monolayers / s (iBiham et al.l 
2001). For these conditions and a grain temperature of T = 18K, 
the results of the KMC simulations are compared to the master 
equation steady-state efficiency in Eq. (112b using either the ap- 
proximation to A given by Eq. © or the conventional A w a/S in 
Fig. [2] Evidently, the agreement is perfect using the approxima- 
tion of A given in Eq. ||6), while the conventional approximation 
is in significant disagreement. 



5. Generalization to multi-species reactions 

Our method has been formul ated for single-species re actions ex- 
clusively. As speculated by lLohmar & Krugl d2006h . we now 
show that it can easily be generalized to multi-species reactions. 
We use X and Y to denote two different species of particles. First, 
we consider the two-particle problem to find the encounter prob- 
ability pxY of an X and a Y particle on a translationally invari- 
ant lattice. So far, we have been concerned with two particles, 
each with a (undirected) hopping rate a and desorption rate W; 
this situation can be mapped onto a single walker moving and 
desorbing at rates 2a and 2W, respectively. Both rates can then 
be rescaled to their old value s, since only their ratio is relevant 
(lLohmar & Kru g 2006, 2009). For two species X and Y, and the 
corresponding rates ax, ay and Wx, Wy, we map to one walker 
moving with rate ax + ay and desorbing at a rate Wx + Wy- The 
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Fig. 2. Recombination efficiency versus grain size, as obtained 
from KMC simulations (circles), from the master equation using 
Eq. © (solid line) and using the conventional approximation 
a/S (dashed line). 

rescaling argument is not applicable for arbitrary rates. However, 
it is evident that we may simply substitute each a in the single- 
species p (Sect. |2]i by ax + ay now, and each W by Wx + Wy. 
This yields the encounter probability pxy of a given pair of X 
and Y particles. 

Second, we derive the relation between the e ncounter proba- 
bility p xy and the sweeping rate Axy, following lLohmar & Krugl 
(2006). We start from the production rate of XY molecules, 



Rxy = Axy(NxNy), 



(13) 



in an analogous way to the single-species situation. This means 
that the reaction reduces the number of XY pairs at a rate Axy 
(compared to a rate 2A in the single-species case). An argument 
fully analogous to the single-species reaction case infers the re- 
lation to the encounter probability pxy of a given XY-pair. The 
encounter probability is given by the ratio of the rate at which 
the reaction takes away a pair to the overall rate at which a (par- 
ticular) pair is removed, by either the reaction or the desorption 
of either constituent, 



Consequently, 



Pxy = 



1XY 



iXY 



A X y + Wx + Wy 



(W x + W y )pxy 



(14) 



(15) 



1 - Pxy 

is the definition for the reaction rate coefficient of X + Y — > XY. 

Since the right-hand side has the same form as for the single- 
species reaction with the substitution W — > Wx + Wy (and hop- 
ping rates do not occur here), the sum of desorption rates can- 
cels with the factor of the pxy asymptotics just as it does for the 
single-species case. Using Eq. ©, we thus obtain the approxi- 
mation 

1 



ax + fly 



IXY 



-71 



In 



i 



for the reaction rate coefficient of X + Y — > 
1.8456047840.... 



(16) 



XY, where c 



expression is numerically accurate, and easy and efficient to 
implement even in complex gas-grain reaction networks. For 
single-species reactions, the expression was employed in the ex- 
act analytical master equation framework, and comparison with 
KMC simulations shows excellent agreement throughout. The 
results were also generalized to reactions of multiple species. 
We strongly advocate the incorporation of the reaction rate co- 
efficients presented here into all models of interstellar gas-grain 
chemistry. This could be easily achieved and, in stark contrast to 
the commonly used approximation, it provides accurate results 
at virtually no additional computational cost. This is an impor- 
tant step towards an accurate modeling of the fascinating chem- 
istry of interstellar clouds. 
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6. Conclusions 

We have used exact results for the sweeping rate in diffusion- 
mediated reactions on dust grain surfaces to provide an ap- 
proximate formula for the crucial reaction rate coefficient. This 



